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We study how quantum and thermal noise affects synchronization of two optomechanical limit- 
cycle oscillators. Classically, in the absence of noise, optomechanical systems tend to synchronize 
either in-phase or anti-phase. Taking into account the fundamental quantum noise, we find a regime 
where fluctuations drive transitions between these classical synchronization states. We investigate 
how this “mixed” synchronization regime emerges from the noiseless system by studying the classical- 
to-quantum crossover and we show how the time scales of the transitions vary with the effective 
noise strength. In addition, we compare the effects of thermal noise to the effects of quantum noise. 


I. INTRODUCTION 

The field of cavity optomechanics deals with systems where light in an optical cavity couples to mechanical motion 
[1]. Recent experimental progress allows to move on from the investigation of a single optomechanical system to several 
coupled optomechanical systems. There are already first experiments that involve a few mechanical and optical modes 
[2-5] , exploiting them for wavelength conversion, phonon lasing and efficient cooling. Such few-mode optomechanical 
setups have been the subject of an increasing number of theoretical proposals, on topics such as efficient state transfer 
[6], two-mode squeezing [7], back-action evading measurements [8], entanglement [9-11] or Landau-Zener dynamics 
[12, 13]. 

Larger arrays may be implemented using a variety of settings, such as coupled disks [14] (Fig. 1(b)) or optome¬ 
chanical crystal structures [15-20] (Fig. 1(d)). Optomechanical arrays have also attracted attention from a theoretical 
point of view. They have been studied in the context of slowing light [21], Dirac physics [22], reservoir engineering 
[23], artificial magnetic fields for photons [24], heat transport [25], and topological phases of sound and light [26]. Fur¬ 
thermore, multi-membrane systems [27-30] were studied theoretically, considering for instance long-range interactions 
and dynamics. 

Most notably, optomechanical arrays provide a platform to study synchronization of mechanical oscillators. This was 
initially pointed out in Ref. [31]. Synchronization is a well known phenomenon in many different branches of science 

[32] and typically arises whenever there are stable limit-cycle oscillations. However, we note that synchronization¬ 
like phenomena have been recently studied also in the context of linear oscillators dissipating into a common bath 

[33] . Optomechanical systems exhibit a Hopf bifurcation and can be optically driven into mechanical limit-cycle 
oscillations [34-37]. These self-oscillations have also been analyzed theoretically in the quantum regime [38-40]. The 
theoretical description of the synchronization dynamics in optomechanical arrays has initially focussed on the classical 
regime [31, 41]. More recent insights into this regime include the pattern formation of the mechanical phase field in 
larger optomechanical arrays [42]. In the quantum regime, it was found [43] that quantum noise can drive a sharp 
nonequilibrium transition towards an unsynchronized state in an extended array, even for optomechanical systems 
with identical frequencies. Further general insights into quantum synchronization were gained in a model system 
of one van-der-Pol oscillator coupled to an external drive [44] or two coupled van-der-Pol oscillators [45], which can 
serve as a rough approximation to an actual optomechanical system. A number of more recent works have explored 
quantum synchronization on the more conceptual level [46], as well as in various other physical systems, such as 
e.g. trapped atoms and ions [47-50], qubits [51-53], and superconducting devices [54]. The relation of synchronization 
and correlations in the quantum-to-classical transition was studied in a system of coupled cavities containing a non¬ 
linearity [55]. Notably, it is still challenging to define a good measure for quantum synchronization [49, 56, 57]. 

Only recently synchronization of two nanomechanical oscillators in the classical regime was demonstrated experi¬ 
mentally. This has been achieved using optomechanical systems with coupled micro-disks [14] (Fig. 1(b)), as well as 
in an experiment involving an optical racetrack cavity coupled to two mechanical oscillators [58] (Fig. 1(c)), and also 
in a setup using nanoelectromechanical systems [59]. In a recent first step towards larger arrays, up to seven optically 
coupled micro-disks were used to demonstrate the expected phase noise reduction due to synchronization [60]. This 
1/N phase noise reduction with the number of coupled systems N is considered to be one of the main prospects of 
synchronized nanomechanical arrays. Indeed, recognized from the very beginning with the synchronization of pendu¬ 
lum clocks [61], synchronization has the potential to improve time-keeping and frequency stability. Examples where 
different types of synchronization have been applied or suggested for application are for frequency stabilization of 
high power lasers by coupling to a more stable, low power laser [62] and for secure communication in connection with 
chaos [63]. For a more complete overview and also the many applications to biology see e.g. [32, 64, 65]. 
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Figure 1: Synchronization of optomechanical oscillators. Schematics of the setup we study (a) and possible experimental 
implementations (b)-(d). (b) Micro-disk oscillators that support optical whispering gallery modes which couple evanescently to 
each other [14, 60], (c) optical racetrack resonator coupled to two nanomechanical oscillators [58], (d) optomechanical crystal 
structure with two optomechanical cells in the vicinity of each other, allowing for optical and mechanical coupling [20]. 


In this work we study the effects of quantum and thermal noise on the synchronization of two optomechanical 
systems. We focus on a bistable synchronization regime that either exists already in the absence of noise, or is 
induced by it. Bistabilities in quantum systems have been investigated before [69-71] and noise-induced bistabilities 
are known in several other systems in biology and chemistry [66-68], as well as in physics [72-74]. We find and discuss 
noise-induced bistable behaviour now in the context of optomechanical (quantum) synchronization. 

This manuscript is organized as follows: We begin with a brief review of classical synchronization in the absence of 
noise, Sec. II. Then, we introduce our model in Sec. Ill, explain our methods in Sec. IV, and state our main results 
in Sec. V. We note that both quantum and thermal noise lead to similar effects, but start out with the investigation 
of quantum noise effects. Therefore, in Sec. VI, we simulate the full quantum behaviour of the system, in contrast 
to the previous investigation presented in Ref. [43]. We find a regime of “mixed” synchronization with two stable 
synchronization states, and we explore its classical-to-quantum crossover in Sec. VII. In Sec. VIII, we give an overview 
of the different synchronization regimes. Finally, in Sec. IX, we discuss the effects of thermal noise as compared to 
quantum noise. This is important to gauge the potential of observing the quantum noise effects discussed here in 
future experiments. 


II. BRIEF REVIEW: CLASSICAL SYNCHRONIZATION OF OPTOMECHANICAL OSCILLATORS 

In this section we briefly review the concepts of classical synchronization of optomechanical systems in the absence 
of noise. This will set the stage for the discussion of our results on synchronization in the presence of thermal and 
quantum noise. 

A widely studied model for synchronization is the Kuramoto model [75, 76], which describes a set of coupled phase 
oscillators. Each phase oscillator has a phase fa and an intrinsic frequency and couples to the other oscillators 
via the phase difference. For two oscillators, the two corresponding phase equations collapse into one equation for the 
relative phase 8(j) = fa — 0i, 


8(j) = (uj2 — coq) — k sm(5fa). (1) 

The two oscillators are synchronized if their respective phase velocities become equal, fa = fa, i.e. 8(j) = const. From 
Eq. (1) one finds the synchronization threshold: The two oscillators are synchronized if the coupling k exceeds the 
natural frequency difference, \k\ > \u>2 — u>i\. Following this condition, oscillators with identical intrinsic frequencies 
uji are always synchronized. There is only a single stable value of 5(f) in the synchronized regime. For k -A Too, this 
value approaches zero, 8(f) 0, whereas 8(f) —>> tt for k — oo. 

Although synchronization appears in a large range of systems which are very different in terms of microscopic 
parameters, their behaviour can often be captured by effective phase equations of the Kuramoto-type [32, 76]. In the 
context of optomechanics, the mechanical motion of a single optomechanical system near the Hopf bifurcation can 
be described with a phase and an amplitude equation [36]. Starting from these equations, an effective Kuramoto- 
type model for coupled optomechanical systems has been derived [31]. This model describes arrays of arbitrary many 
optomechanical cells with arbitrary intrinsic frequencies. For two oscillators, the phase equations of this model reduce, 
again, to a single equation for the phase difference 
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Figure 2: Threshold of self-sustained oscillations and synchronization in the absence of noise, (a) shows the mechanical energy 
(related to the limit-cycle amplitude) of an optomechanical system as a function of the laser drive and the detuning A. 
The red line indicates the threshold of self-oscillations, T + T op t = 0. The green line separates a region of optical multi¬ 
stability, the dashed white line indicates a region where the self-oscillations show strong amplitude modulation, (b) and (c) 
show schematic pictures of the classical, noiseless synchronization regimes of two coupled optomechanical systems. In (b) we 
use a small mechanical coupling K/Qi = 0.05 and different mechanical resonance frequencies, 5Q = 0.075 to give an overview 
about all possible synchronization regimes. In (c) we use K/Q± = 0.15 and 5Q = 0 instead, which significantly changes the 
phase diagram. The dashed lines mark regions with even more complex behaviour, not necessarily showing synchronization. 
[Parameters: mechanical damping T/fli = 0.015, optical damping k/Qi = 0.3, optomechanical coupling go/Qi = 0.3] 


8(j) = (uj 2 — kff) — 2Si sin(£</>) — 452 sin(2 Sc/)). (2) 

Here, Si and S 2 are effective parameters depending on the microscopic parameters of the underlying optomechanical 
systems [31, 43]. Note that the 5 ' 1 -term was added to the model only later [43], and accounts for the change of 
the intrinsic frequencies Ui with the mechanical oscillation amplitude. Recently, the resulting full Hopf-Kuramoto 
model has been used to study pattern formation in 2D arrays of optomechanical systems [42]. In contrast to the 
original Kuramoto model Eq. (1), the optomechanical Hopf-Kuramoto model Eq. (2) includes a term that involves 
sin(2£</>). Rewriting Eq. (2) in terms of an effective potential, 5<j) = —U'(5(f)) [31], this term corresponds to the 
appearance of a second minimum close to S(f) = tt which can co-exist with the minimum close to Sc/) = 0. This 
allows optomechanical systems to synchronize not only in-phase (0-synchronization), 5(f) -T 0, but also anti-phase 
(7r-synchronization), 8(j) 7 r. For the effective potential there are three different possibilities, depending on the 
parameters S\ and 52: (i) It has a single minimum close to Sc/) = 0 which leads to 0-synchronization only, (ii) it has 
a single minimum close to Scj) = tt which leads to 7r-synchronization only, or (iii) both minima appear simultaneously 
in the effective potential, such that the initial conditions determine whether 0- or 7r-synchronization occurs. These 
three regimes are schematically shown in Fig. 2(b) and (c) and discussed in Sec. III. 

In general, for optomechanical arrays an effective potential for the phases & does not exist [42]. However, in the 
case of two oscillators only, the system can be described with one degree of freedom, i.e. the relative phase Sq i>, and an 
effective potential for 5(f) can always be constructed. Below, we make use of the existence of this effective potential to 
give an intuitive understanding of synchronization even in the presence of noise. 


III. MODEL 

We now introduce the model that we investigate throughout the rest of this manuscript. We study two optome¬ 
chanical systems which are mechanically coupled, see Fig. 1(a). Our goal is to analyze the synchronization behaviour 
in the presence of quantum and thermal noise. Each optomechanical system consists of a driven optical mode (dj) 
coupled to a mechanical mode (bj) via radiation pressure. In a frame rotating at the laser drive frequency ujl , the 
Hamiltonian of each optomechanical system is 
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Hj = -hAa^dj + - hgoOjdj (bj + + fia L (a] + a^). (3) 

Here, A = ujl —uj c denotes the detuning from the cavity resonance gj c , ft is the resonance frequency of the mechanical 
mode, go denotes the optomechanical single photon coupling strength, and cul is the laser driving strength. Note that 
the optical and the mechanical systems experience damping at a rate k and T, respectively, and these are described 
by adding to the Hamiltonian a system-bath coupling in the usual manner: = Hj + H^ dlss \ 

Driving a single optomechanical system with a blue-detuned laser causes self-oscillations of the mechanical resonator 
when the optomechanically induced negative damping r opt is larger than the intrinsic damping of the oscillator T 
[36]. In Fig. 2(a) we show this threshold of self-oscillations. Throughout this work we only consider parameters such 
that the single optomechanical systems are above this threshold. Self-oscillations at A < 0 occur due to the static 
optomechanical shift which leads to an effective blue detuning, i.e. A e ff > 0. These limit-cycle oscillations, in the 
absence of noise, can effectively be described by a fixed amplitude and a phase and are treated as a prerequisite 
for synchronization throughout this work. We consider two self-oscillating optomechanical systems that are coupled 
mechanically with strength K, such that the total Hamiltonian of the system (except for the dissipative part) reads 


H tot = ^2 Hj- hK (&! + 61) (i>2 + 61) . (4) 

3 = 1,2 

Experimentally, this coupling between the mechanical oscillators can also be mediated by an optical coupling, 
cf. Fig. l(b)-(d). In recent experiments [14, 58, 60], a single joint optical mode was employed to couple the me¬ 
chanical oscillators. However, this mode served a dual purpose in creating the limit cycles via a blue-detuned drive 
and providing the coupling. Using an additional, independently driven optical mode for the coupling would allow to 
tune the effective (optically induced) coupling independently from the laser drive used to create the limit cycles. 

In experiments, the typical mode of operation is to have the two optomechanical oscillators at slightly different 
intrinsic mechanical frequencies. These start out un-synchronized but can synchronize upon changing some parameter 
(e.g. the laser drive strength, the detuning, or potentially the coupling). This is schematically shown in Fig. 2(b), 
where we indicate the synchronization regimes in the absence of noise. In accordance to the Hopf-Kuramoto model, 
cf. Sec. II, there are unsynchronized regions and three different synchronization regimes: (i) O-synchronization, (ii) 
7r-synchronization, and (iii) classical bistable synchronization where the type of synchronization depends on the initial 
conditions. Note that for different intrinsic mechanical frequencies, 5Q = D 2 — Cli 7 ^ 0, the relative phase 5(f) is not 
exactly 0 or tt but only close to one of these values and varies within the synchronization regime. 

However, in the presence of noise it is already interesting to investigate the behaviour even for identical frequencies. 
In particular, for large-scale optomechanical arrays of identical oscillators, it has been found that there is a synchro¬ 
nization transition as a function of noise strength [43]. Moreover, in the present article we will focus on noise-induced 
transitions between various synchronization states. The observation of this physics does not rely on whether there is 
an actual synchronization transition at lower values of the coupling. Therefore, in most of our analysis, we will focus 
on identical systems, i.e. we assume all the parameters to be equal in both systems. In Fig. 2(c) we schematically 
show the synchronization regimes for identical optomechanical systems and for a larger mechanical coupling K than 
in Fig. 2(b), but still in the absence of noise. It is important to note that both the mechanical detuning 5Q and the 
coupling K have an influence on this diagram: Not synchronized regions can become synchronized and synchronization 
regime borders are shifted. 

We will comment on the dynamics of two coupled optomechanical oscillators with different frequencies in Sec. VIII. 


IV. METHODS 

In this work we use Langevin equations and quantum jump trajectories to study the system described by Hamil¬ 
tonian (4) in the presence of (quantum) noise. Here we want to briefly present both approaches and discuss their 
respective advantages and problems. 

Most of our results are computed with semi-classical Langevin equations. They are obtained by first deriving 
quantum Langevin equations from Hamiltonian (4) using input-output theory [ 88 ]. We then adopt a semi-classical 
approach by turning the quantum Langevin equations into classical Langevin equations for the complex amplitudes 
Oij and /3j, where the noise terms mimic the quantum-mechanical zero-point fluctuations. This can be understood as 
a variant of the “truncated Wigner approximation”. The semi-classical equations are: 
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di = (iA - |) ai + (/?i + PI) ~ iot L - V^aim, 

/ p\ (5) 

$1 = ~ f ^ ^ ) Pi + ^0 |<M1 2 + iK@2 — VrPlin- 

The corresponding equations for the second optomechanical system can be obtained from Eqs. (5) by exchanging 
the indices 1 «—» 2. Here, aqi n (t) and represent the optical and mechanical input noise, given by Gaussian 

stochastic processes. Since complex numbers commute, they obviously cannot correctly fulfill the input-output quan¬ 
tum noise correlators, (a^ in (t)dji n (t')) = 0 and (aji n (t)aj in (t / )) = S(t — t'). Instead, aj- in and /3j[ n are made to mimic 
quantum noise by fulfilling (aji n (t)a* in (t')} = = 5(t — t ')/2 (and likewise for (3j[ n for T = 0). This 

approach allows to study also large parameter ranges with relatively low computational effort. 

Deep in the quantum regime it is initially not clear that these semi-classical Langevin equations describe the correct 
physical behaviour. In order to verify the qualitative effects observed with Langevin equations, we also present a few 
results obtained with quantum jump trajectories [77, 78], i.e. an “unraveling” of the Lindblad master equation. 
Applying this method, the fully quantum system is simulated on an appropriately truncated Hilbert space. Notably, 
it allows to work with wave functions, in contrast to the formalism of the Lindblad master equation which requires 
the use of density matrices. Hence, simulations of larger Hilbert spaces become feasible. In addition, quantum jump 
trajectories give access to additional observables, for instance the full counting statistics. For these reasons, quantum 
jump trajectories have been applied in many different contexts. In the field of quantum synchronization, such methods 
have been used to study synchronization of qubits [51]. In the field of cavity optomechanics, they have been employed 
for instance to discuss QND measurements, photon statistics and single photon optomechanics [79-82]. Another 
motivation is the recent experimental detection of individual phonons in optomechanical systems [83] - a step towards 
monitoring full quantum jump trajectories in experiments. 

To explain this approach [77, 78], let us consider for a moment photon decay in cavity 1. The “unraveling” discussed 
here corresponds to the physical setup of placing a single photon detector at the output port of the cavity. At each time 
step 5t, the probability of a single photon to leak out of the cavity (at temperature T = 0) is given by p = K>5t(d\di). In 
the case of a photon loss it is detected at the output port and the wave function is updated |T(£ + St)) = d\ |T(£)). If 
no photon was lost to the environment, the wave function evolves in time according to a (non-Hermitian) Hamiltonian 
that is obtained by adding a term —ih{n/2)d\d\. This additional term accounts for the information gained about the 
system by not observing a photon [84]. In both cases the state |T(£ + St)) has to be normalized before proceeding to 
the next time step. The treatment of photon loss in cavity 2 and of phonon losses works analogously. This simulation 
approach naturally accounts for quantum fluctuations. 

Although it would be favourable to use quantum jump trajectories throughout the whole study, this approach is only 
computationally feasible whenever the needed Hilbert space remains sufficiently small. This leads to severe limitations 
in the choice of parameters and especially gives no access to the full classical-to-quantum crossover that is studied 
in Sec. VII. In Ref. [38] it has been shown for a single optomechanical system that semi-classical Langevin equations 
produce good agreement with the full quantum theory. A systematic comparison for two coupled optomechanical 
systems is not possible, since the number of required Hilbert space dimensions is squared as compared to the single 
optomechanical system. We comment on this in Sec. VIIA (cf. Fig. 4(a) and the discussion below). 

Note that quantum jump trajectories serve here as a numerical approach only. In experiments, single photon and 
phonon detection is not necessary. Instead, the mechanical oscillators could be measured using standard homodyning 
techniques. 


V. MAIN RESULTS 

In this section we briefly state our main results, which will be discussed in the following sections. Most notably, 
it is known even from the classical theory that two coupled optomechanical oscillators can be in either one of two 
synchronization states (with a phase difference near 0 or near 7r). We find a regime of “mixed” synchronization, where 
transitions between 0- and 7r-synchronization occur (Sec. VI). These transitions are driven by (quantum or thermal) 
noise and cannot be found in the classical, noiseless situation. The average residence times in the two synchronization 
states can differ and their ratio varies with the system parameters (Sec. VII). Investigating the classical-to-quantum 
transition, we find that mixed synchronization can evolve from two different regimes in the classical, noiseless limit: 
(i) there are already two stable synchronization states but in the absence of noise there are no transitions, (ii) there 
is only one stable synchronization state and only the presence of noise leads to a second stable solution. Although 
the first sections are devoted to the investigation of quantum noise effects, we note that we find similar effects for 
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Figure 3: Multistability in optomechanical quantum synchronization. Distribution of the relative phase 5(f) in different syn¬ 
chronization regimes, a typical sample of a corresponding quantum jump trajectory starting in the steady state, and a sketch 
of the corresponding effective potential, (a) shows O-synchronization, (b) shows 7r-synchronization, and (c) shows mixed syn¬ 
chronization. A rotating wave approximation for the mechanical coupling has been used. [Parameters: (a) mechanical coupling 
K/Q = 0.3, mechanical damping r/D = 0.015; (b) K/Q = 0.15 , T/D = 0.01; (c) K/Vt — 0.15, T/Q = 0.015; other parameters 
are: optical damping k/Q = 0.3, laser driving strength /D = 0.3, optomechanical coupling go/n — 1, and optical detuning 
A/D - 0.15. ] 


thermal noise acting on the mechanical resonator. However, quantitative differences remain due to the different nature 
of the noise source (Sec. IX). We find that quantum noise effects should dominate over thermal noise effects if the 
optomechanical cooperativity is sufficiently large, and a large value of go is not necessarily required. 


VI. MULTISTABLE QUANTUM SYNCHRONIZATION 

First, we start by analyzing two coupled optomechanical systems, cf. Fig. 1(a), deep in the quantum regime. 
Quantum jump trajectories are used to initially investigate the full quantum dynamics. In the following, we consider 
a small-amplitude limit cycle and a large single-photon coupling strength, go/ft = 1. This ensures that quantum 
fluctuations can potentially have a large impact on the system’s dynamics. Furthermore, small photon and phonon 
numbers are necessary to keep the numerical simulations tractable, since they determine the size of the truncated 
Hilbert space. 

To study synchronization we focus on the relative phase 5(f) = 02 — (f>i between the optomechanical oscillators. 
Classically, 5(f) allows to identify synchronization ( 5(f) = const) and distinguish the different synchronization regimes (0- 
and 7r-synchronization). To extract the relative phase at each time step of a quantum jump trajectory, we use that the 
motion of the uncoupled self-oscillators in the absence of noise can effectively be described by ( bj)(t ) ~ 
with an initial random phase (f)Qj. Thus, the correlator ~ e~ l6< ^ can be interpreted as a measure for the 

relative phase 5(f). In Fig. 3 we show the probability density of the relative phase obtained from the distribution of 
5(f) and a typical sample of the corresponding phase trajectory. We find a 0-synchronized regime, Fig. 3(a), where 
the relative phase is predominantly close to 5(f) = 0. Similarly, we find a 7r-synchronized regime, Fig. 3(b), for 
different parameters. As expected, noise prevents perfect synchronization and thus the corresponding maximum in 
the probability density has a finite width. Likewise, the trajectories show fluctuations around either 5(f) = 0 or 5(f) = tt. 
These two synchronization regimes have an analog in the classical, noiseless limit, cf. Sec. II. In addition, we find 
another regime where the probability density of the phase difference has maxima close to both 5(f) = 0 and 5(f) = 7r, 
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see Fig. 3(c). The corresponding trajectory shows that in this regime transitions between the two synchronization 
states occur. We call this regime “mixed” synchronization. This is in contrast to the classical, noiseless result, where 
the system does not change its synchronization state with time, not even in the regime of bistable synchronization. 

Similar to the classical, noiseless system [31] we can understand these results in terms of an effective potential for 
the relative phase 5(f). This is illustrated in the bottom row of Fig. 3. It offers an intuitive understanding of our 
results. The noise makes the phase fluctuate around the stable point(s) near 5(j) = 0 (and)or 5(j) = i r. In the case of 
mixed synchronization the effective potential has two minima, near 5(f) = 0 and Scj) = 7r, and quantum noise drives 
transitions between those two states. The probability to be in either the 0- or 7r-synchronized state is given by the 
area of the corresponding maximum in the probability density. It is associated to the depth of the minimum in the 
effective phase potential. The ratio between the probabilities of the two states can assume arbitrary values, depending 
on the parameters of the system. Note that the absence of two distinct peaks in the probability density does not 
necessarily mean that the system explores the region around one synchronization state only. In fact, sometimes the 
trajectories themselves can reveal short stretches of phase dynamics in the vicinity of the other state. Nevertheless, 
if the fraction of time spent in that other state is short, a second peak will not be visible. 

From this interpretation in terms of an effective potential, it is clear that the classical bistable synchronization 
regime where two potential minima already exist (but no transitions can occur), turns into a mixed synchronization 
regime (showing transitions) in the presence of noise. However, our analysis in the following section reveals that mixed 
synchronization can also appear for parameters where classically there is only a single stable synchronization state. 


VII. CLASSICAL-TO-QUANTUM CROSSOVER 

At the moment, the single-photon coupling strength go is still comparatively small in almost all experiments. As 
a consequence, quantum effects have only been observed in the linearized regime where only Gaussian states are 
produced. In the most promising cases [85, 86], go/& can take values up to 10 -2 and go/Q up to above 10 -4 . 
Much larger values have been reported for experiments with cold atoms (up to around go/ft ~ 1) [87], but these 
do not operate in the “good cavity limit”, i.e. one has ft Q in those experiments, precluding the observation of 
single-photon strong coupling effects. Nevertheless, as experiments are approaching the single-photon strong coupling 
regime, they will gradually see increasingly strong effects of quantum fluctuations even in non-linear dynamics. In 
this section, it is our aim to explore the crossover between the classical regime (small go /ft) and the quantum regime 
(large go /ft) with respect to optomechanical quantum synchronization. We will focus our investigations mostly on the 
mixed synchronization regime which is the most interesting one, as we can have noise-induced transitions between 
the synchronization states. In this section, we will disregard thermal noise, i.e. we assume temperature T = 0, such 
that only quantum noise is present. The “classical” regime we are discussing here is therefore the noiseless limit of 
the classical equations of motion. We will later remark on the effects of thermal noise (Sec. IX). 

To explore the classical-to-quantum crossover, we want to effectively vary h while making sure to keep all the 
classical predictions unchanged. Notably, the optomechanical coupling strength go = duj c /dx depends on h. For the 
simplest case of a Fabry-Perot cavity of length L with one static and one movable mirror the optomechanical coupling 
is g 0 = u c (0)xzpf/L \/h, where xzpf = ^/S/2 mfl denotes the zero-point fluctuations of the mechanical oscillator of 
mass m. As discussed in Ref. [38], the “quantum parameter” go/ft ^ Vh can thus be varied to effectively change the 
quantum noise strength. This implies that all classical (^-independent) parameters (ft/f2, T/Q, A/fi, K/Q, go a L /^ 2 ) 
are kept fixed while go is modified. To see that the quantum parameter has indeed this anticipated effect, it is very 
helpful to rescale the amplitudes ay and /3j, Eqs. (5), such that they tend to a well-defined finite value in the classical 
limit #o/ft —>■ 0 [38]. This can be achieved by defining ay = gootj and /3j = go/3j. While |aq| 2 gives the number of 
photons (a “quantum-mechanical” quantity), the rescaled version |aq| 2 = g^ |aq| 2 ~ h |aq| 2 is proportional to the 
energy hui c |aq| 2 inside the cavity (i.e. a classical quantity). A similar argument applies to ftj. With this rescaling, we 
find the following equations: 


bq = (iA - aq + idq (fti + j - ia L - V^iin, 

^ / p\ _ _ (6) 

Pi = — ( i£l + — ) Pt + i |bq | iKfo — VT/3u n . 

Here dz, = goo^L is the rescaled laser-driving amplitude that we keep fixed while varying go . The important observation 
here is that go has been completely eliminated from the equations and now only appears in the strength of the quantum 
noise: we now have (dji n (£)d* in (£')) = g^{t — t')/2 (and likewise for /3j[ n ) which indeed vanishes in the classical limit 

of go Vh —y 0. 




quantum parameter go /k 



Figure 4: Quantum-to-classical crossover, (a) Synchronization measure as a function of the quantum parameter go/a, using 
Langevin equations (black line), classical (noiseless) Langevin equations (black circle, go/n = 0), and with quantum jumps (red 
triangles), (b) and (c) show trajectories of the relative phase in the quantum regime go/n = 1 (blue square) and for go /ft = 0.4 
(green dot). Other parameters are as in Fig. 3(c), a rotating wave approximation has been used for the mechanical coupling. 


If we consider mechanical thermal noise at finite temperature T (as discussed in Sec. IX), we have (/3ji n (t)/3* in (t')) = 
#o( n th + 1/2 )S(t — t'). Here, n t h denotes the thermal occupancy of the bath coupled to the mechanical oscillator. The 
product #o n th ~ ftrith becomes independent of H in the classical limit >> HQ, where Hn t h ~ ksT/Q and is 
Boltzmann’s constant. In summary, the rescaled equations (6) nicely show how an increase of go can indeed be viewed 
solely as an increase of the strength of “quantum noise” in our system. 

Note that if ^ 0 we have to increase the laser driving strength oll to keep oll = go&L constant which means that 
the total light energy circulating inside the optical cavity is constant. Due to H 0, this corresponds to an increasing 
number of photons inside the cavity and thus increases drastically the size of the Hilbert space necessary for the full 
quantum simulation. Therefore, quantum jump trajectories are not suitable to explore the full quantum-to-classical 
crossover and we have to apply Langevin equations. 


A. Synchronization as a Function of Quantum Noise Strength 

In all studies of synchronization, one needs to select suitable quantities that measure the degree of synchronization. 
In this context, it is important to note that for any finite noisy system (subject to quantum and/or thermal noise), 
there is no sharp synchronization transition and correspondingly no unambiguous measure that displays nonanalytic 
behaviour at any parameter value. Before proceeding to our results, we summarize and discuss the synchronization 
measures adopted here which have to be combined to obtain a full picture: (i) the probability density of 8(f), (ii) the 
average phase factor and (iii) individual trajectories. 

The probability density of the relative phase 5(j) is either mostly flat (no synchronization) or, as shown in Fig. 3, 
centered predominantly around 0 or 7r, or it may have two peaks, according to the synchronization regime. In the 
following, we aim to compress the information contained in the phase distribution into one quantity and calculate the 

normalized correlator C = $$ 2 )/\J (^ 2 ^ 2 ) ~ Its real value, Re[C] = (cos£</>), distinguishes the three 

different synchronization regimes: (i) (cos 5<fi) ~ 1 for O-synchronization, (ii) (cos 5(p) ~ —1 for 7r-synchronization, (iii) 
intermediate values of (cos 5<p) for mixed synchronization. However, this measure has its limitations: When Scj) is more 
or less evenly distributed (no synchronization) (cos 5(j)) ~ 0, this cannot be distinguished from a mixed synchronization 
situation where almost equal time is spent in the 0- and 7r-synchronized states. Furthermore, even in the absence 
of synchronization (and even in the noiseless case) the phase 5(j) may spend an increased amount of time around 
certain values. This leads to a finite value of (cos £</>), and similarly would also show up in the phase distribution. A 
solution to this problem is to simultaneously look at a part of the corresponding trajectory, where synchronization can 
easily be distinguished from an unsynchronized state. Instead, one could start to use more complicated correlators, 
e.g. ( b\b\b 2 b 2 ) ~ e - l26( t >. This correlator allows to distinguish unsynchronized states from synchronized states, but an 
additional measure is needed to distinguish 0- from 7r-synchronization. Note that the imaginary part of the above 
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Figure 5: Residence times of synchronization states, (a) and (b) show the distribution of residence times in the 0-synchronized 
and 7r-synchronized state. The yellow histogram shows the distribution for go/n = 0.4, the green histogram gives the result 
deeper in the quantum regime, go/ft = 1, for comparison. Note the different scale of the time axis in (a) and (b). The 
corresponding typical time scale To or T n is obtained from a fit and is shown in (c) as a function of the quantum parameter 
go/ft. Results are obtained using Langevin equations, with the following parameters: optical damping k,/Q = 0.3, mechanical 
damping T/Q = 0.015, detuning A/Q = —1/30, mechanical coupling K/Q = 0.15, gootL/Q 2 = 0.09 = const, and a rotating 
wave approximation for the mechanical coupling term. 


defined correlator, Im[C] = (sin Sc/)), can be used as well. However, in the special case of identical optomechanical 
systems lm[C] ~ 0 due to the symmetry of the system. 

In Fig. 4(a) we show how (cos 5(f)) varies as a function of the quantum parameter go / k. We chose parameters that lead 
to a mixed synchronization regime for larger values of go/ft. In the deep quantum regime (go/ k 1) the probability 
Po to find the system in the 0-synchronized state is larger than the probability P\ to find the phase around 7r, such 
that (cos 5(f)) > 0. Going towards smaller values of go/ ft, the ratio P^/Po increases, such that eventually (cos 5(j)) < 0. 
Finally, we should reach the classical (noiseless) limit, when go /k —> 0. It turns out that, for the parameters adopted 
here, the classical solution always ends up in the 7r-synchronized state, independent of initial conditions. This implies 
that there is only one minimum in the effective potential. We conclude that the system has turned from a mixed 
synchronization regime into a purely 7r-synchronized regime as the quantum parameter was reduced. This cannot be 
understood in the simple picture of a noise-independent phase potential. We will discuss this kind of behaviour in 
more detail later on (Sec. VII C). 

Note that the calculations for Fig. 4 have been performed using Langevin equations; although in the deep quantum 
regime, two data points were also acquired with quantum jump trajectories (red triangles). They are shifted as 
compared to the Langevin results, but show the same trend. We expect that the difference between the Langevin and 
quantum jump results decreases for smaller quantum parameters go/ k, as it is the case for a single optomechanical 
system [38]. Since smaller go/ft require a significantly larger Hilbert space for the quantum jump simulations, we 
cannot compute this for our coupled system. For large quantum parameter go/ ft ^ 1, qualitative differences between 
the full quantum model and the Langevin equations have been already observed in Ref. [38] as well. Especially a shift 
of the detuning A was reported, that could be determined numerically also for our system. Taking this detuning shift 
into account would improve the agreement of our results, although differences remain. Here, we show the uncorrected 
outcomes of both approaches. 


B. Residence Times in the Mixed Synchronization Regime 

The measure (cos 5(j)} quantifies the fraction of total time spent in 0-synchronized parts as compared to 7r- 
synchronized parts. However, it does not provide any information about the rate of transitions between the two 
synchronization states. Based on the effective potential picture, one might expect the transition rates to be deter¬ 
mined by the barrier height and the noise strength. In particular, for larger effective noise strengths go/ k, we expect 
more frequent transitions. This behaviour is qualitatively visible in Fig. 4(b) and (c). However, as concluded in the 
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relative phase 5(f) 


Figure 6: Appearance of bistability. The distribution of the relative phase as a function of the quantum parameter go/n. 
Parameters are as in Fig. 5, but with mechanical coupling K/Q = 0.1. 


previous section, the potential picture is not sufficient to explain all observations. Thus, we now turn to a quanti¬ 
tative analysis and discuss how the transition rates between 0- and 7r-synchronized states (i.e. the typical residence 
times r) change during the classical-to-quantum crossover. We extract the fluctuating residence times from the phase 
trajectories and obtain their distribution. The results are shown in Fig. 5(a) and (b) for the 0- and 7r-synchronized 
states, for two different quantum parameters go/ft. In all cases the probability densities decay exponentially with 
time, ^ e -r / r , and the average residence time f is obtained from a fit to the distribution. The extracted average 
residence times To and for the two states are shown in Fig. 5(c) as a function of the quantum parameter. Note 
that the ratio of residence times equals the ratio of probabilities, ro/r n = Po/P Nevertheless, the dependence of the 
times To, r n on go /k reveals new information. 

We have chosen parameters such that at go /k = 1 both 0-synchronized and ^-synchronized parts have almost equal 
average residence times. This corresponds to (cos 6(f)) « 0 and Po ~ P^. Furthermore, in the classical limit go/ft = 0 
the system is ^-synchronized only. When the classical limit go /k -a 0 is approached, we find that both To and r n 
increase. As expected, r n increases much faster than To and eventually diverges for go/ft -A 0, as the system gets 
trapped forever in the 7r-state. In contrast, To increases first when decreasing go/ft, but then saturates at a finite 
level. Such a behaviour is unexpected based on the simple phase potential picture, where a fixed potential would 
imply diverging residence times for both states in the noiseless limit. The behaviour observed here hinges on the fact 
that the synchronization regime switches from “mixed” to “t” as one reduces the quantum parameter (i.e. reduces 
the quantum noise). The observations would change significantly for different parameters, where the system always 
stays in the mixed regime, for any value go /k > 0. Then, one expects the simple picture of a fixed phase potential to 
be approximately correct and both residence times to diverge as the noise is becoming weaker. 


C. Noise-induced Synchronization Bistability 

In the previous sections we have explained that some basic features of the classical-to-quantum crossover, like the 
increase of the residence times with decreasing quantum parameter, can be understood as effects of a decreasing 
quantum noise strength. The decrease of noise strength leads to less frequent transitions across an energy barrier 
in the effective potential. However, we made also less easily explained observations: (i) the reverse in the order 
of To and r n as go is being reduced, (ii) the saturation of To at low noise levels, (iii) the disappearance of stable 
0-synchronization (for the applied parameters) in the classical limit go /k — 0. Whereas (i) and (ii) could originate 
from more complicated potential shapes with a combination of broad and narrow minima, (iii) suggests that the 
effective potential itself changes when the quantum parameter is varied. In Fig. 6 we show how the distribution 
of 5(f) evolves as a function of go/ft. For very small values of go /«, there is only a single peak close to 5(f) ~ 7r, in 
accordance with the single stable solution of the classical limit. While increasing the quantum parameter, this peak 
is first broadened. The increasing quantum noise strength allows the system to explore more of the effective potential 
around the minimum. A significant accumulation close to 5(f) ~ 0 appears only for rather large quantum parameters, 
signaling the appearance of a second stable solution, i.e. a second minimum in the effective potential at 5(f) ~ 0. 

In addition to the above described appearance of a second stable solution, there are also parameter regions where 
already in the classical regime the effective potential has minima at both 5(f) = 0 and 5(f) = it (classical bistable 
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synchronization). In this case quantum noise naturally drives transitions between the two synchronization states as 
soon as it is added to the description of the system. The number of observed transitions then naturally depends on 
both the noise strength as well as the potential shape. 


VIII. OVERVIEW OF SYNCHRONIZATION REGIMES 

In the previous sections we have shown examples of the different synchronization regimes in the presence of quantum 
noise and studied the properties of mixed synchronization in more detail. In the following, we map out the different 
synchronization regimes as a function of the mechanical coupling K and the quantum parameter go/K. Furthermore, 
we also discuss the case of detuned mechanical oscillators. 

Figure 7(a) shows the synchronization measure (cos Sc/)} for resonant oscillators. We have indicated the synchro¬ 
nization regimes (which are not sharply delineated). For the applied parameters, we find O-synchronization for large 
mechanical coupling K in both the quantum and classical regime. In the classical, noiseless limits (cos Sc/)) 1, 

indicating less fluctuations around the synchronization state. In contrast, at smaller mechanical coupling, we find 
more complicated behaviour: there is mixed synchronization for go/K ~ 1 , while the classical limit go /k —> 0 selects 
either 0— or tt— synchronization, depending on the mechanical coupling K. The “pixelated” region in Fig. 7(a) indi¬ 
cates that the system is multistable even in the classical limit. There, the residence times have become so large that 
the system is stuck in a random synchronization state depending on initial conditions and the transient behaviour. 
Notably, the closer the system is to a border of synchronization regimes in the classical limit (this can be seen in 
Fig. 7(a) for small go/K when K is varied), the smaller the noise strength (go/K) that is needed to lead to mixed 
synchronization. As an example, in the ”middle” of the classically tt- synchronized regime (at about K/YL\ « 0.1) 
similar mixed synchronization as compared to the system close to the regime border (K/Yt\ < 0.15) appears only 
for larger values of the quantum parameter. An exception is of course the classically bistable regime, where mixed 
synchronization appears naturally as soon as there is noise. An interesting feature appears close to K/Yl ^ 0.18, 
where the measure (cos Sc/)) shows a sharp dip in the middle of a 0 -synchronized region. We suspect a non-linear 
resonance, since the oscillator trajectory Xj(t) is no longer simply sinusoidal and period-doubling is observed. At the 
same time, the oscillation amplitude increases. For larger values of the mechanical coupling A, the trajectories are 
simply sinusoidal again, with the same frequency as for coupling strengths below the feature. 

Up to now, we have studied the ideal case of identical mechanical oscillators. We now turn to the case where the 
mechanical oscillators have slightly different resonance frequencies, i.e. = YI 2 — fii 7 ^ 0. This is a typical situation 
in experiments, since fabrication inaccuracies lead to deviations between two systems. Figure 7(b) shows how c)U 
affects synchronization. A finite SQ corresponds to a tilted effective potential in the classical, noiseless limit, where 
a finite threshold for synchronization appears [31]. This classical threshold is indicated in Fig. 7(b) with a dashed 
line. Similar behaviour is visible for quantum parameters go /k > 0. However, for large go /k it is not possible to 
determine the onset of synchronization, because the measure (cos S(/)) cannot distinguish between no synchronization 
and mixed synchronization. Instead, we now also have to analyze the trajectories in more detail, which reveal mixed 
synchronization for sufficiently large mechanical coupling. A significant deviation from the classical threshold cannot 
be observed at the given resolution. We expect the threshold to slightly increase for larger quantum parameter go/K. 
At the same time, however, the threshold is also smeared out due to quantum noise. 

Here, we have chosen to show the dependence of the synchronization regimes on the mechanical coupling strength 
K and the quantum parameter go/K. Note that other parameters influence the synchronization type as well. In Fig. 2 
we have already seen that the laser driving strength oil and the detuning A influence the classical synchronization 
regimes. Also the mechanical damping Y can affect the observed synchronization, cf. Fig. 3b and c. However, first 
of all it already influences the limit cycles of individual optomechanical systems by modifying the threshold to self- 
sustained oscillations. In addition, V has an influence on the mechanical noise strength. Note that, when changing 
these parameters, care has to be taken to remain on a stable limit cycle for each optomechanical system. 


IX. THERMAL NOISE 

So far, we assumed zero temperature environments for both the optical and the mechanical mode. In this section 
we investigate the effects of thermal mechanical noise on synchronization. 

For our study we use the Langevin equations (5) with modified mechanical noise correlators to account for the 
coupling of the mechanical oscillators to a finite temperature bath, (/3j{ n (tfPjinit')) = (PjinitfPjinit')) = ( n th + 
1 / 2 )S(t — £ 7 ), where n t h denotes the thermal occupancy of the bath. Hence, both quantum and thermal noise are 
included. For optical frequencies in the visible spectrum the effective thermal occupation of the optical bath is very 


12 


oscillators resonant 


oscillators detuned 


oscillators resonant 



0.5 1.0 

quantum parameter go/n 


0.5 

quantum parameter 


0 


-1 


50000 


thermal noise Tith 


o 

o 

U1 


Figure 7: Synchronization regimes versus mechanical coupling and noise. The measure (cos Sc/)) as a function of the quantum 
parameter go/n (i.e. effective quantum noise strength) and of the mechanical coupling strength K/Q i, for identical optome¬ 
chanical systems (a) and for systems with mechanical frequency detuning SQ/Q i = (Q 2 — Hi)/Qi = 0.05 (b). In (c), we display 
the effects of thermal noise nth- In that case, we show (cos 5<j)) as a function of coupling K and environmental thermal phonon 
number n t h, with a constant go /ft = 0.01 and for zero S£l. The grey region in (a) and (b) was not simulated. The dashed line 
in (b) is the synchronization threshold as observed in the classical, noiseless limit go/ft = 0. Parameters are as in Fig. 5. 


small. Thus, the assumption of an optical bath at zero temperature is valid and the optical input noise terms are not 
modified. 

In Fig. 7(c) we show an overview of the synchronization regimes as a function of thermal noise n t h- Here, we chose 
a comparatively small quantum parameter go/hi = 0.01 in order to observe mainly effects due to thermal noise. For 
small n t h the results are similar to Fig. 7(a) at small go/hi. In both cases the influence of the quantum or thermal noise 
is still weak. For increasing thermal noise strength n t h, we find qualitatively the same behaviour as for increasing 
quantum noise. However, quantitative differences appear with increasing n t h- Even though a mixed synchronization 
regime can develop for both quantum and thermal noise, the evolution of the relative weight of both synchronization 
states is different. 

In the following, we want to estimate the critical thermal noise strength n* h at which thermal and quantum noise 
should have a comparable effect. At lower temperatures, quantum noise will dominate. The main source of quantum 
fluctuations is the laser shot noise (for the parameters explored here). Thus, we estimate the effect of optical quantum 
fluctuations on the mechanical oscillator, with the (symmetrized) shot-noise spectrum evaluated at the mechanical 
resonance frequency £1 [1, 89], 




1 
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( %o 
V^ZPF 



K 

(k/2) 2 + (^ + A) 2 


+ ( K /2)2 + (ft- A) 2 ) • 


(7) 


We expect similar effects from quantum and thermal noise if the shot-noise spectrum at the mechanical resonance fre¬ 
quency £1 becomes equal to the thermal force spectrum, Spp(Q) = S l FF . The thermal noise spectrum at temperatures 


ksT Ml is 


Sp F = 2mTk B T, (8) 

where T is the temperature of the thermal mechanical bath. Setting S)^ F and S l FF equal, we find (in the resolved 
sideband regime 9 at A ^ 9): 


n* h = k B T/Ml = C/2, (9) 

where we used the optomechanical cooperativity [1] 

C = (10) 

This approach suggests that observing quantum noise phenomena does not necessarily require a large go/hi and very 
low temperatures. Instead, if the cooperativity is sufficiently large (comparable to values that enable ground-state 
cooling, C > n t h)j quantum noise should dominate the behaviour of the system even in the presence of thermal noise. 
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However, depending on parameters, we find large deviations from this simple expectation. In these cases, the real 
shot noise spectrum is no longer well described by the weak-coupling expression of S'f^(O) given above, and the 
actual noise strength may have a much larger value. Consequently, the transition between behaviour dominated by 
quantum noise vs. that dominated by thermal noise takes place at much larger values of n* h than those predicted by 
Eq. (9). In other words, it should be even easier to observe quantum noise in optomechanical synchronization than 
the naive ansatz would lead one to expect. 

Here, we don’t observe that quantum noise can be exactly mapped to thermal noise. This is already evident in the 
rescaled Langevin equations (6). Physically, the thermal force spectrum acting on the mechanical oscillator is flat in 
frequency, whereas the optical shot-noise spectrum is frequency dependent and is also modified by the dynamics of 
the system. 


X. CONCLUSIONS 

We have investigated the effects of quantum and thermal noise on two coupled optomechanical limit-cycle oscillators. 
One usually expects that noise prevents strict synchronization, i.e. exact phase locking and a sharp transition to syn¬ 
chronization. Here we have shown that fluctuations additionally drive transitions between 0- and n- synchronization, 
i.e. the two synchronization states that can appear in the absence of noise. We have discussed the residence times of 
these states and observed a smooth crossover between different synchronization regimes. Finally, we have compared 
the effects of quantum and thermal noise. We have argued that it should be possible to experimentally reach the 
regime where quantum noise dominates. This should happen when the optomechanical cooperativity is large enough 
for ground state cooling. 

For further investigations it would be useful to identify a measure that can genuinely distinguish between an 
unsynchronized regime and 0-, i r- and mixed synchronization. Finally, it will be very interesting to extend the 
insights obtained here to large optomechanical arrays. 
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